% % % % % % % % % % % % % % % % % % % % % % % % % % 
%  task.m - first version written on Aug 5, 2011  %
% % % % % % % % % % % % % % % % % % % % % % % % % % 

% First, start log of output
diary('log.out');

% Load some useful variables from previous search results.
load('workspace.mat');

% Set parameters for the search. Start with initial parameter values, then set chosen model specification and grid size parameters. 

init        = Optparams;
init(9:10)  = [0.0003 0.0003];
span        = [0.020  0.100  0.0015  0.0015  0.0005  0.075  0.050  0.050  0.000  0.000  0.000  0.050  0.050];
nK     		= 100;
nL     		= 100;
nA     		= 2;
nF     		= 2;
extraA 		= 10;
Lmode  		= 'LaborNet';
FC     		= 'Fixed';
Spec        = 1;

fprintf('----------------------------------------- \n')
fprintf('  Initialize parameter search \n')
fprintf('----------------------------------------- \n')
fprintf('nK     = %3.0f             \n',nK)
fprintf('nL     = %3.0f             \n',nL)
fprintf('nA     = %3.0f             \n',nA)
fprintf('nF     = %3.0f             \n',nF)
fprintf('extraA = %3.0f             \n',extraA)
fprintf('----------------------------------------- \n')
fprintf('Lmode  = LaborNet        \n')
fprintf('FC     = Fixed             \n')
fprintf('Search over full set of parameters. \n') 
fprintf('----------------------------------------- \n \n')


if strcmp(Lmode,'LaborGross') == 1
   momlist   = [10 14 18 23 25 90 91 92 94 69 89];
elseif strcmp(Lmode,'LaborNet') == 1
   momlist   = [10 14 18 23 25 90 91 93 95 69 89]; 
end

% Start the annealing search with given inputs:
[Criterion Optparams Optmom results NSE] = annealing(init, span, nK, nL, nA, nF, extraA, Lmode, FC, momlist, Spec);
load('empmom.mat','EY')
fprintf('------------------------------------ \n')
fprintf('     Results: Search completed   \n')
fprintf('------------------------------------- \n')
fprintf('Estimates:    Parameter     St.err. \n')
fprintf('------------------------------------- \n')
fprintf('bK             %s       %s  \n', results{1,1}, results{1,2})
fprintf('bL             %s       %s  \n', results{2,1}, results{2,2})
fprintf('alphaK         %s       %s  \n', results{3,1}, results{3,2})
fprintf('alphaL         %s       %s  \n', results{4,1}, results{4,2})
fprintf('alphaKL        %s       %s  \n', results{5,1}, results{5,2})
fprintf('RSloss         %s       %s  \n', results{6,1}, results{6,2})
fprintf('rho            %s       %s  \n', results{7,1}, results{7,2})
fprintf('sigma          %s       %s  \n', results{8,1}, results{8,2})
fprintf('sigdK          %s       %s  \n', results{9,1}, results{9,2})
fprintf('sigdL          %s       %s  \n', results{10,1}, results{10,2})
fprintf('ro             %s       %s  \n', results{11,1}, results{11,2})
fprintf('meK            %s       %s  \n', results{12,1}, results{12,2})
fprintf('meL            %s       %s  \n', results{13,1}, results{13,2})
fprintf('------------------------------------- \n \n')
fprintf('------------------------------------ \n')
fprintf('Criterion = %6.3f \n', Criterion)
fprintf('------------------------------------ \n \n')

fprintf('--------------------------------------------------------- \n')
fprintf('     Moments           Empirical   Simulation  Diff\n')
fprintf('--------------------------------------------------------- \n');
fprintf('Corr(IK(t),IK(t-2))      %7.4f     %7.4f     %7.4f  \n', EY(10), Optmom(10), abs(EY(10)-Optmom(10)))
fprintf('Corr(HL(t),HL(t-2))      %7.4f     %7.4f     %7.4f  \n', EY(14), Optmom(14), abs(EY(14)-Optmom(14)))
fprintf('Corr(Ygr(t),Ygr(t-2))    %7.4f     %7.4f     %7.4f  \n', EY(18), Optmom(18), abs(EY(18)-Optmom(18)))
fprintf('--------------------------------------------------------- \n')
fprintf('Std(IK(t))               %7.4f     %7.4f     %7.4f  \n', EY(23), Optmom(23), abs(EY(23)-Optmom(23)))
fprintf('Std(HL(t))               %7.4f     %7.4f     %7.4f  \n', EY(25), Optmom(25), abs(EY(25)-Optmom(25)))
fprintf('--------------------------------------------------------- \n')
fprintf('p(I/K<0)                 %7.4f     %7.4f     %7.4f  \n', EY(90), Optmom(90), abs(EY(90)-Optmom(90)))
fprintf('p(I/K=0)                 %7.4f     %7.4f     %7.4f  \n', EY(91), Optmom(91), abs(EY(91)-Optmom(91)))
fprintf('p(H/L=0)                 %7.4f     %7.4f     %7.4f  \n', EY(92), Optmom(92), abs(EY(92)-Optmom(92)))
fprintf('p(dL=0)                  %7.4f     %7.4f     %7.4f  \n', EY(93), Optmom(93), abs(EY(93)-Optmom(93)))
fprintf('p(I/K=0;H/L=0)           %7.4f     %7.4f     %7.4f  \n', EY(94), Optmom(94), abs(EY(94)-Optmom(94)))
fprintf('p(I/K=0;dL=0)            %7.4f     %7.4f     %7.4f  \n', EY(95), Optmom(95), abs(EY(95)-Optmom(95)))
fprintf('--------------------------------------------------------- \n')
fprintf('Skewness(I/K)            %7.4f     %7.4f     %7.4f  \n', EY(69), Optmom(69), abs(EY(69)-Optmom(69)))
fprintf('Skewness(H/L)            %7.4f     %7.4f     %7.4f  \n', EY(89), Optmom(89), abs(EY(89)-Optmom(89)))
fprintf('--------------------------------------------------------- \n \n')


[criterion simmom IK_use HL_use K_use L_use Y_use V_use A_use Optk Optl lKmin lLmin lKmax lLmax dK_use dL_use] = simulation(Optparams, nK, nL, nA, nF, extraA, Lmode, FC, 'Search',momlist);

fprintf('---------------------------------------------------- \n')
fprintf(' Simulated impact of adjustment costs (yearly avg.) \n')
fprintf('---------------------------------------------------- \n')
fprintf(' Sales value                       = %7.4f \n', mean(mean(mean(Y_use))))
fprintf(' Wage bill                         = %7.4f \n', 1*mean(mean(mean(L_use))))
fprintf(' Capital replacement               = %7.4f \n', 1*.1*mean(mean(mean(K_use))))
fprintf(' Convex capital adjustment costs   = %7.4f \n', (1/2)*Optparams(1)*mean(mean(mean(IK_use)))^2*mean(mean(mean(K_use))))
fprintf(' Convex labor adjustment costs     = %7.4f \n', (1/2)*Optparams(2)*mean(mean(mean(HL_use)))^2*mean(mean(mean(L_use))))
fprintf(' Fixed capital adjustment costs    = %7.4f \n', Optparams(3)*mean(mean(mean(abs(IK_use)>0.015))))
fprintf(' Fixed labor adjustment costs      = %7.4f \n', Optparams(4)*mean(mean(mean(abs(HL_use-dL_use)>0.015))))
fprintf(' Interrelated adjustment costs     = %7.4f \n', Optparams(5)*mean(mean(mean((abs(IK_use)>0.015).*(abs(HL_use-dL_use)>0.015)))))
fprintf(' Resale loss on capital sold       = %7.4f \n', Optparams(6)*mean(mean(mean((IK_use<0).*IK_use.*K_use)))*(-1))
fprintf('---------------------------------------------------- \n')

diary('off')
save('workspace.mat')

exit
